# Gilardi, Fabrizio, "The Temporary Importance of Role Models for Women's Political Representation", American Journal of Political Science
# Code to replicate Figure SI5 (Effect of an additional woman elected in the same municipality on the number of women candidates)
# gilardi@ipz.uzh.ch, 2014-06-23

# Set working directory
setwd("../Data/")

# Load data
d <- read.csv("Frauenvertretung-kantonale-Parlamente-2.csv", stringsAsFactors=F)
e <- read.csv("einwohner-1970-2009.csv", stringsAsFactors=F)


n.k <- nrow(d)
n.yrs <- ncol(d)-1
yrs <- c(seq(1971, 2011, 4), 2013)

dd <- as.data.frame(matrix(NA, nrow=n.k*n.yrs, ncol=4))

dd[,1] <- rep(d$Kanton, n.yrs)

j <- 1
for(i in 1:n.yrs){
	dd[j:(j+n.k-1),2] <- rep(yrs[i], n.k)
	j <- j + n.k
}

j <- 1
for(i in 1:n.yrs){
	dd[j:(j+n.k-1),3] <- d[,(i+1)]
	j <- j + n.k
}


yrs.e <- c(1,3,7,11,15,19,23,27,31,35,39,41)

e <- e[,yrs.e]

j <- 1
for(i in 1:(n.yrs-1)){
	dd[j:(j+n.k-1),4] <- e[,(i+1)]
	j <- j + n.k
}

colnames(dd) <- c("canton", "year", "pcent.w", "pop")

de <- unique(dd$canton)[c(1:9,11:20)]

dd$german <- 0
dd$german[(dd$canton %in% de) == T] <- 1

dd$pcent.w[dd$year < 1972 & (dd$canton=="St. Gallen" | dd$canton=="Uri" | dd$canton=="Schwyz" | dd$canton=="Graubünden" | dd$canton=="Nidwalden" | dd$canton=="Obwalden")] <- NA

dd$pcent.w[dd$year < 1989 & (dd$canton=="Appenzell A. Rh.")] <- NA
dd$pcent.w[dd$year < 1990 & (dd$canton=="Appenzell I. Rh.")] <- NA

dd$yr.w <- NA

dd$yr.w[dd$canton==unique(dd$canton)[1]] <- 1970
dd$yr.w[dd$canton==unique(dd$canton)[2]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[3]] <- 1971 
dd$yr.w[dd$canton==unique(dd$canton)[4]] <- 1972
dd$yr.w[dd$canton==unique(dd$canton)[5]] <- 1972
dd$yr.w[dd$canton==unique(dd$canton)[6]] <- 1972
dd$yr.w[dd$canton==unique(dd$canton)[7]] <- 1972
dd$yr.w[dd$canton==unique(dd$canton)[8]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[9]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[10]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[11]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[12]] <- 1966
dd$yr.w[dd$canton==unique(dd$canton)[13]] <- 1968
dd$yr.w[dd$canton==unique(dd$canton)[14]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[15]] <- 1989
dd$yr.w[dd$canton==unique(dd$canton)[16]] <- 1990
dd$yr.w[dd$canton==unique(dd$canton)[17]] <- 1972
dd$yr.w[dd$canton==unique(dd$canton)[18]] <- 1972
dd$yr.w[dd$canton==unique(dd$canton)[19]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[20]] <- 1971
dd$yr.w[dd$canton==unique(dd$canton)[21]] <- 1969
dd$yr.w[dd$canton==unique(dd$canton)[22]] <- 1959
dd$yr.w[dd$canton==unique(dd$canton)[23]] <- 1970
dd$yr.w[dd$canton==unique(dd$canton)[24]] <- 1959
dd$yr.w[dd$canton==unique(dd$canton)[25]] <- 1960
dd$yr.w[dd$canton==unique(dd$canton)[26]] <- 1977

dd$yrs.since.w <- dd$year - dd$yr.w


dd$cntn <- NA

dd$cntn[dd$canton==unique(dd$canton)[1]] <- "ZH"
dd$cntn[dd$canton==unique(dd$canton)[2]] <- "BE"
dd$cntn[dd$canton==unique(dd$canton)[3]] <- "LU" 
dd$cntn[dd$canton==unique(dd$canton)[4]] <- "UR"
dd$cntn[dd$canton==unique(dd$canton)[5]] <- "SZ"
dd$cntn[dd$canton==unique(dd$canton)[6]] <- "OW"
dd$cntn[dd$canton==unique(dd$canton)[7]] <- "NW"
dd$cntn[dd$canton==unique(dd$canton)[8]] <- "GL"
dd$cntn[dd$canton==unique(dd$canton)[9]] <- "ZG"
dd$cntn[dd$canton==unique(dd$canton)[10]] <- "FR"
dd$cntn[dd$canton==unique(dd$canton)[11]] <- "SO"
dd$cntn[dd$canton==unique(dd$canton)[12]] <- "BS"
dd$cntn[dd$canton==unique(dd$canton)[13]] <- "BL"
dd$cntn[dd$canton==unique(dd$canton)[14]] <- "SH"
dd$cntn[dd$canton==unique(dd$canton)[15]] <- "AR"
dd$cntn[dd$canton==unique(dd$canton)[16]] <- "AI"
dd$cntn[dd$canton==unique(dd$canton)[17]] <- "SG"
dd$cntn[dd$canton==unique(dd$canton)[18]] <- "GR"
dd$cntn[dd$canton==unique(dd$canton)[19]] <- "AG"
dd$cntn[dd$canton==unique(dd$canton)[20]] <- "TG"
dd$cntn[dd$canton==unique(dd$canton)[21]] <- "TI"
dd$cntn[dd$canton==unique(dd$canton)[22]] <- "VD"
dd$cntn[dd$canton==unique(dd$canton)[23]] <- "VS"
dd$cntn[dd$canton==unique(dd$canton)[24]] <- "NE"
dd$cntn[dd$canton==unique(dd$canton)[25]] <- "GE"
dd$cntn[dd$canton==unique(dd$canton)[26]] <- "JU"


yrs[10]

i <- 10
out <- lm(pcent.w ~ log(pop) + yrs.since.w, data=dd, subset=c(year==yrs[i]))
summary(out)

can <- dd$cntn[dd$year==yrs[i] & is.na(dd$pcent.w)==F]
act <- dd$pcent.w[dd$year==yrs[i] & is.na(dd$pcent.w)==F]


# Figure for 2011

par(mar=c(3.6,3.5,2,1), mgp=c(2.5,0.8,0), cex.axis=1.2, cex.lab=1.2, font.main=1)
plot(predict(out), act, type="n", xlab="Predicted % of women in parliament", ylab="Actual % of women in parliament", axes=F)
axis(1)
axis(2)
text(predict(out)[-1], act[-1], can[-1], col="gray", cex=1.2)
text(predict(out)[1], act[1], can[1], cex=1.2)
abline(0,1)
#title(yrs[i])
legend("topleft", legend=yrs[i], cex=1.5, bty="n")


# Figures for each year

for(i in 1:11){

f <- paste("cantons-", i, ".pdf", sep="")	
pdf(file=f, paper="special", width=6, height=6)

out <- lm(pcent.w ~ log(pop) + german, data=dd, subset=c(year==yrs[i]))
summary(out)

can <- dd$cntn[dd$year==yrs[i] & is.na(dd$pcent.w)==F]
act <- dd$pcent.w[dd$year==yrs[i] & is.na(dd$pcent.w)==F]

par(mar=c(3.6,3.5,2,1), mgp=c(2.5,0.8,0), cex.axis=1.5, cex.lab=1.5, font.main=1)
#plot(predict(out), act, type="n", xlab="Predicted % of women in parliament", ylab="Actual % of women in parliament", axes=F)
plot(predict(out), act, type="n", xlab="", ylab="", axes=F)
axis(1)
axis(2)
text(predict(out)[-1], act[-1], can[-1], col="gray", cex=1.5)
text(predict(out)[1], act[1], can[1], cex=1.5)
abline(0,1)
#title(yrs[i])
legend("topleft", legend=yrs[i], cex=1.5, bty="n")

dev.off()

}


